3x3 Householder QR Demo

This demo constructs a $3\times 3$ QR factorization using Householder reflectors.

In [10]:
import numpy as np

import numpy.linalg as la
In [11]:
n = 3



e1 = np.array([1,0,0])

e2 = np.array([0,1,0])

e3 = np.array([0,0,1])



A = np.random.randn(n, n)

A
Out[11]:
array([[ 1.06427586,  1.03785412,  1.19643376],
       [-1.06405176,  0.48693156, -0.83335032],
       [-0.20227028, -1.10857722,  0.27986689]])

Householder reflector:

$$I-2\frac{vv^T}{v^Tv}$$

Choose $v=a-\|a\|e_1$.

In [12]:
#clear

a = A[:, 0]

v = a-la.norm(a)*e1



H1 = np.eye(3) - 2*np.outer(v, v)/(v@v)
In [13]:
A1 = H1 @ A

A1
Out[13]:
array([[  1.51848691e+00,   5.33870203e-01,   1.38523070e+00],
       [  3.40401588e-16,  -6.93719959e-01,  -3.91067578e-01],
       [  9.99443798e-17,  -1.33301245e+00,   3.63942360e-01]])

NB: Never build full Householder matrices in actual code! (Why? How?)

In [14]:
#clear

a = A1[:, 1].copy()

a[0] = 0

v = a-la.norm(a)*e2



H2 = np.eye(3) - 2*np.outer(v, v)/(v@v)
In [15]:
R = H2 @ A1

R
Out[15]:
array([[  1.51848691e+00,   5.33870203e-01,   1.38523070e+00],
       [ -2.45801147e-16,   1.50272073e+00,  -1.42307423e-01],
       [ -2.55820086e-16,   1.61523465e-16,   5.14914060e-01]])
In [16]:
Q = np.dot(H2, H1).T

la.norm(np.dot(Q, R) - A)
Out[16]:
8.5458331069102901e-16
In [ ]: